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Abstract The Tokay Gecko, Gekko gecko (Linnaeus, 1758) is widely distributed in Asia and there have been concerns 
regarding locally decreasing populations due to overexploitation for traditional Chinese medicine. Previous studies of 
the genetic relationships of G. gecko populations included few populations from Thailand. Here we investigated the 
phylogeographic patterns of G. gecko from different regions in Thailand using mitochondrial cytochrome b sequences. 
Phylogenetic analyses revealed two lineages: one (Lineage A) comprising populations from Laos, Vietnam, and 
Thailand; and a second (Lineage B) comprising three genetically distinct groups within Thailand alone. Some Thai 
populations were found to have both lineages represented within them. Highly significant genetic differentiation (Fsr) 
showed geographic population structuring in Lineage B, indicating limited gene flow among groups in Thailand. 
Although G. gecko has a wide distribution and is well adapted to human habitation, the observed genetic structure could 
potentially be explained by geographic barriers such as mountain ranges. In Lineage A, our study provided primary 
phylogeographic evidence for lineage mixture that might be a result of human-mediated transport. Future research 
should include more extensive sampling across the geographic distribution of G. gecko and a landscape genetics 
approach could be applied for conservation planning. 
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1. Introduction et al., 2012). Cryptic species are two or more distinct 
species erroneously classified (and hidden) under a single 
Southeast Asia has a unique geological history that has species name (Bickford et al., 2007). In some cases 
contributed to high regional biodiversity (Sodhi et al., two or more members of a cryptic species complex can 
2004; Woodruff, 2010). To understand the diversification 


of organisms and historical biogeography of this 


occur syntopically, making identification challenging and 
confounding the interpretation of results from studies 
based on specimens from these populations (McLeod, 
2010; Stuart et al., 2006). Moreover, the genetic structure 
of populations can be affected by human activities such 


region, studies on adaptation and evolution of reptiles 
in response to the geographical and climatic changes 
are of great interest (e.g., Brown et al., 2012; Heinicke 
et al., 2012; Siler et al., 2013), particularly studies 
exploring cryptic diversity of reptiles at both species 


as habitat fragmentation and human-mediated transport 
(Santos et al., 2018; Templeton et al., 2001). In these 
instances, genetic divergence can be used to determine 


end population: leveli tee Brown ar ah, AULA ema whether individuals and populations belong to the same 


: or separate evolutionary lineages (Ferguson, 2002) and 
Corresponding author: Dr. Anchalee Aowphol, from Department of 


Zoology, Faculty of Science, Kasetsart University, Bangkok, Thailand, 
with her research focusing on systematics and ecology of amphibians 
and reptiles in Southeast Asia. 

E-mail: fsciacl@ku.ac.th 


Received: 30 November 2018 Accepted: 1 July 2019 


can indicate the level of connectivity among populations 
(e.g., Blair et al., 2013; Dubey et al., 2011). 

Among vertebrates, reptiles are good models for 
phylogeographic studies because of their global 
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distribution (Pincheira-Donoso et al., 2013; Uetz, 2000), 
and a high degree of genetic divergence (Lin et al., 2010). 
Very few phylogeographic studies of reptiles have been 
conducted in Thailand. Lukoschek et al. (2011) and 
Saijuntha et al. (2017) used mitochondrial DNA (mtDNA) 
to study the phylogeography of the Mekong Mud Snake 
(Enhydris subtaeniatus) and Blue Crested Lizard (Calotes 
mystaceus), respectively. The results of these studies have 
suggested that genetic structure of reptile populations in 
Thailand might be affected by geographic barriers such as 
mountain ranges and river drainage systems. 

The Tokay Gecko Gekko gecko (Linnaeus, 1758) 
belongs to the family Gekkonidae and has a wide 
distribution in Asia (Rösler, 2005; Smith, 1935; Uetz 
et al., 2018). It has been introduced to many regions of 
the world and is considered an invasive species in the 
United States and Brazil (Meshaka et al., 1997; Rocha, 
2015). The Tokay Gecko has long been heavily exploited 
for Chinese traditional medicine, and this has led to 
localized population declines (Bauer, 2009) in parts of 
Bangladesh, Indonesia, and Thailand (Caillabet, 2013). 
The taxonomic status of G. gecko remains controversial 
and it is reported to comprise two morphotypes (black- 
spotted and red-spotted geckos). These two morphotypes 
occur allopatrically and have high genetic divergence 
based on mitochondrial and nuclear DNA (Qin et al., 
2007, 2012; Wang et al., 2012, 2013), and different 
chromosome characteristics (Qin et al., 2012). Rösler 
et al. (2011) revalidated the taxonomic status of the 
black-spotted gecko as Gekko reevesii (Gray, 1831) 
based on morphological, molecular, and zoogeographic 
evidence. In addition, the two morphotypes have different 
advertisement calls (Yu et al., 2010) and ecological 
data showed niche differentiation (Zhang et al., 2014), 
corroborating Rösler et al. (2011)’s revalidation of 
G. reevesii in China and northern Vietnam. There 
are currently two recognized subspecies based on 
morphology, namely G. gecko gecko (Linnaeus, 1758), 
and G. gecko azhari Mertens, 1955. Gekko g. gecko is 
widely distributed, ranging from northeast India, southern 
China, and Southeast Asia, whereas G. g. azhari is found 
only in Bangladesh (Rösler et al., 2011; Smith, 1935; 
Rösler, 2001, 2005). This study considers only samples 
from outside of Bangladesh and therefore we refer the 
study taxon (G. g. gecko) herein as G. gecko. 

Although G. gecko has a wide distribution range in 
Asia, phylogeographic analyses of this species from some 
geographical areas, including Thailand, have been limited 
(Kongbuntad et al., 2016; Qin et al., 2007, 2012; Wang et 
al., 2012, 2013; Zhang et al., 2006). This study examines 
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the genetic structure of G. gecko populations in Thailand. 
We used mitochondrial cytochrome b (cyt b) DNA 
sequence data to assess the mitochondrial DNA diversity 
and phylogeographic patterns of G. gecko populations in 
Thailand, and compare Thai samples to specimens from 
Laos, Vietnam, and southern China. 


2. Materials and Methods 


2.1. Sample collection Field work was conducted 
from 2010-2018. A total of 114 samples were collected 
from 17 populations in the residential areas throughout 
Thailand where G. gecko occurs as a human commensal 
(Table 1; Figure 1). Tail tips of G. gecko were clipped, 
and geckos were returned to where they were captured. 
Liver tissues were collected from individuals that were 
kept as voucher specimens. All tail tip and liver tissues 
were preserved in 95% ethanol and stored at —20 °C for 
genetic analyses. Voucher specimens were fixed in 10% 
formalin and subsequently transferred to 70% ethanol. 
All tissues and voucher specimens were deposited in the 
herpetological collection, Zoological Museum, Kasetsart 
University (ZMKU). 


2.2. DNA extraction, PCR amplification, and 
sequencing Total genomic DNA was extracted from 
tail tip or liver tissues preserved in 95 % ethanol using 
a DNeasy Blood and Tissue Kit (Qiagen, Hilden, 
Germany) according to manufacturer’s protocol. Part 
of the mitochondrial cyt b gene was amplified by 
polymerase chain reaction (PCR) using the combination 
of forward and reverse primers from Wang et al. (2012) 
(Table 2). Amplification conditions included initial DNA 
denaturation at 94°C for 5 min; followed by 35 cycles of 
94°C for 1 min, 50°C annealing temperature for 50 s, and 
72°C for 1 min; with a final extension at 72°C for 10 min. 
PCR products were purified using a QIAquick purification 
kit (Qiagen, Hilden, Germany) and sequenced in both 
forward and reverse directions using the PCR primers. 
The sequencing reactions were analyzed on a 3730 DNA 
Analyzer (Applied Biosystems, CA, USA) by Macrogen 
Inc. (Seoul, Korea). Sequences were checked by eye, 
and aligned using MUSCLE implemented in Geneious 
R11 (Biomatters, Ltd., Auckland, New Zealand). All 
DNA sequences obtained in this study were deposited in 
GenBank (Accession numbers MK598399-MK 598441). 


2.3 Data analyses 

Phylogenetic reconstruction 

Additional cyt b sequences of G. gecko from Phonsavan, 
Laos (GenBank: EF174461, EF640157), Langson, 
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Table 1 Locality, samples sizes (n), Number of haplotype, MtDNA lineage, and haplotype distribution of Gekko gecko used in the present 


study. 
No. (Province, Distr) Code Latitude Longitude Th) type (NR). lineage (No individua) o 
1 Loei, Phu Ruea LE 17.422684 101.459448 7 2 B2 H20 (1), H21 (6) 
2 Khon Kaen, Waeng Noi KK 15.810563 102.415894 7 1 B2 H19(7) 
3 plas Ratchasima, Wang Nam NR 14.509586 101.931036 19 6 A2, B2 ee Hot) 10.8); AM 
4 Chanthaburi, Kang Hang Meaw CT 13.135393 101.97300 4 3 A2 _H4 (1), H5 (2), H6 (1) 
5 Chonburi, Sri Racha CB 13.168862 101.087946 1 1 B3 H28(1) 
6 Saraburi, Sao Hai SR 14.550146 100.847831 17 8 B2 E (1), H14 
7 Phitsanulok, Nakhon Thai PL 17.020493 100.925163 2 2 A1,A2 H2 (1), H5 (1) 
8 Uthai Thani, Lan Sak UT 15.578954 99.413981 4 4 B3 H23 (1), H24 (1), H25 (1), H26 (1) 
9 Kanchanaburi, Thong Pha Phum KB 14.667833 98.598191 7 4 B3 H32 (4), H41 (1), H42 (1), H43 (1) 
10 Ratchaburi 1, Mueang RB1 13.580977 99.831455 8 5 Al, B3 a (USS) H34 (1), H36 (4), H39 
11 Ratchaburi 2, Suan Phueng RB2 13.478048 99.252775 3 3 B2, B3 H10 (1), H11 (1), H39 (1) 
12 Phetchaburi 1, Cha Am PB1 12.885258 100.015093 8 3 B3 H37 (1), H38 (1), H40 (6) 
13 Phetchaburi 2, Kaeng Krachan PB2 12.804437 99.574896 3 3 B3 H29 (1), H30 (1), H35 (1) 
14 aa? Khiri Khan, Bang PK 11.213513 99.513924 16 3 B3 H31 (14), H34 (1), H36 (1) 
15 Ranong, Mueang RN 9.853019 98.621854 2 1 B3 H27(2) 
16 Surat Thani, Ban Ta Khun ST 8.955667 98.805028 2 1 Bl H22 (2) 
17 Songkhla, Hat Yai SK 6.940250 100.254389 1 1 A2 H3(1) 


Vietnam (GenBank: EF640155, EF640156), and Guangxi, 
China (GenBank: EF174462) were downloaded from 
GenBank and included in the phylogenetic analyses. 
Gekko reveesii from Guangxi and Yunnan, China 
(GenBank: EF640150-EF640154, EF640179) and 
Gekko vittatus Houttuyn, 1782 (GenBank: NC008772) 
were chosen as outgroups. Phylogenetic relationships 
among haplotypes were estimated using both maximum 
likelihood (ML) and Bayesian inference (BI). A ML tree 
was constructed in the IQ-TREE web server (Trifinopoulos 
et al., 2016) using the Bayesian information criterion 
(BIC). The HKY+F+G4 model was selected as the best-fit 
model of evolution for all codon positions. A BI analysis 
was performed using the CIPRES Science Gateway portal 
V. 3.3 (Miller et al., 2015) with default priors. The best- 
fit substitution model (GTR+I+G) was selected using the 
Akaike Information Criterion (AIC) as implemented in 
jModelTest 2.1.10 (Daribba et al., 2012). Phylogenetic 
reconstruction was performed running Metropolis- 
coupled Markov chain Monte Carlo sampling with four 
chains for 10 million generations and trees were sampled 


every | 000 generations. The first 25% of trees were 
discarded as burn-in and the 50% majority-rule consensus 
of the remaining trees was constructed to calculate the 
posterior probabilities of nodes. Stationarity was checked 
in Tracer v1.7.1 (Rambaut et al., 2018) to ensure that 
effective sample sizes were above 200 for all parameters. 


Genetic diversity and population genetic analyses 

Haplotypes were extracted using DnaSP v6.12.06 (Rozas 
et al., 2017). Population genetic diversity among and 
within sampled populations i.e., haplotype diversity 
(h; Nei, 1987), nucleotide diversity (z; Nei and Tajima, 
1981), the numbers of haplotypes (Nh), and the numbers 
of unique haplotypes (Vy) were calculated using DnaSP 
v6.12.06. Uncorrected pairwise sequence divergences 
(p-distances) between paired phyloclades were calculated 
in MEGA v7 (Kumar et al., 2018). We constructed a 
median joining network using Network v4.6.0 (Fluxus 
Technology Ltd., England) to determine the haplotype 
relationships and calculated the genetic differentiation 
among gecko populations as pairwise fixation indices 
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Table 2 Primers for PCR amplification and sequencing from Wang 
et al. (2012). 


Primer Sequences (5’—3’) 

Forward CACCAAAACCAGTAGTCCGA 
Forward GACCTTCCAACACCATCAAAC 
Reverse CAAGGCCAGTGTATTTGATGT 
Reverse GGGACAAGTAATGGGCACTT 
Reverse GAGCCCCATTTCTGGTTTAC 


(Fst) in Arlequin v3.11 (Excoffier et al., 2005). Genetic 
structure among and within geographic regions was 


assessed using the Analysis of Molecular Variance 


98° E 101°E 
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(AMOVA) approach implemented in Arlequin. 


Population demographic history 

The possibility of demographic expansion of G. 
gecko populations in Thailand was tested based on the 
mismatch distribution of pairwise differences (Rogers 
and Harpending, 1992) in Arlequin. The Harpending’s 
raggedness index (r) and sum of squared deviations 
(SSD) were calculated to compare observed and expected 
distributions under the expansion model. We calculated 
a neutrality test using Fu’s F,(Fu, 1997) and Tajima’s 
D (Tajima, 1989) in Arlequin with 1 000 permutations. 
Significant negative values of Fu’s F,and Tajima’s D 
would indicate demographic expansion. 
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Figure 1 Map showing sampling localities in Thailand. The details of each locality are shown in Table 1. Blue dots indicate samples of 
Gekko gecko in Thailand. Red dots represent samples of G. gecko in Laos, Vietnam, and China. Black dots represent samples of G. reevesii in 


China. 
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Figure 2 Phylogenetic tree, representing relationship among haplotypes and closely related species. Values above or below branch 
correspond to maximum likelihood bootstrap support and Bayesian posterior probability. 
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3. Results 


3.1. Genetic diversity An alignment of | 125 base pairs 
(bp) of cyt b sequences was obtained. There were 201 
polymorphic sites, including 186 parsimony informative 
sites and 15 singleton variable sites. We identified 43 
haplotypes within 114 individuals from 17 localities in 
Thailand. The haplotype diversity (h) and nucleotide 
diversity (z) of each locality are presented in Table 3. The 
overall haplotype diversity was high (0.960 + 0.007), but 
nucleotide diversity was low (0.04076 + 0.00253). 


3.2 Phylogenetic analyses and haplotype relationships 

Phylogenetic analyses (BI and ML) of haplotypes showed 
similar tree topologies. Gekko gecko populations, formed 
a monophyletic group sister to G. reevesii . The BI tree 
(Figure 2) shows that G. gecko is separated into two major 
clades, A and B. Two subclades were recovered within 
clade A. Subclade A1 contained two individuals from PL 
and RB1, whereas subclade A2 consisted of individuals 
from NR, CT, PL, and SK, together with populations 
from Laos (Phonsavan), Vietnam (Langson), and China 
(Guangxi). Clade B comprised three subclades: B1 (ST 


Table 3 Genetic diversity of each phylogenetic clade. 
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population from southern Thailand), B2 (populations 
from northeastern and central (SR) Thailand, and two 
individuals from RB2), and B3 (populations from western 
and southern Thailand, and one individual from CB). 
However, the relationship among these three subclades 
was only moderately supported (BI posterior probability 
= 0.70, ML bootstrap = 72). 

The median-joining network showed intraspecific 
genetic structure, comprising two major lineages (Clades 
A and B) within populations of G. gecko in Thailand 
(Figure 3) that were consistent with the phylogenetic tree 
topology. Most haplotypes were unique to one individual 
(60.47%). Five haplotypes (H5, H10, H11, H34, and 
H36) were shared among populations. Clade A was 
highly divergent from clade B and contained individuals 
from different regions in Thailand. Most individuals in 
clade A carried a unique haplotype, but HS was shared by 
individuals from CT, NR, and PL. In subclade B2, H10 
was the most abundant haplotype and occurred in geckos 
from three localities (NR, SB, and RB2), while H31 was 
the most abundant haplotype in subclade B3 and was 
unique to individuals from the PK population. 


Subclade Sample Size No. of haplotypes No. of unique haplotype Haplotype diversity Nucleotide diversity 
(n) (N) (Nu) (h+ SD) (t+ SD) 

Al 2 2 2 1.000 + 0.500 0.00178 + 0.00089 
A2 7 4 3 0.714 40.181 0.00279 + 0.00136 
Bl 2 1 1 Z — 

B2 51 15 7 0.907 + 0.018 0.00675 + 0.00059 
B3 52 21 14 0.899 + 0.028 0.00971 + 0.00097 
All 114 43 26 0.960 + 0.007 0.04076 + 0.00253 


Table 4 Uncorrected pairwise genetic divergence (%) between phylogenetic clades for mitochondrial cytochrome b dataset. 


Clade Gekko reevesii Al Bl B2 B3 
(n= 6) (n=2) (n=7) (n= 2) (n=51) (n= 52) 
n 0.25 
Gekko reevesii (0.09-0.53) 
Al 8.56 0.18 
(8.36-8.80) (0.18-0.18) 
A2 8.99 2.62 0.28 
(8.71-9.33) (2.58-2.76) (0.00-0.89) 
8.59 8.18 8.36 
Bl (8.44-8.80) (8.18-8.18) (8.27-8.44) ee 
B2 10.59 8.80 9.22 5.45 0.67 
(10.13-11.29) (8.44-9.33) (8.80-10.04) (5.07-6.04) (0.00-1.33) 
B3 10.18 9.10 9.31 3.66 5.44 0.97 
(9.60-10.93) (8.80-9.60) (8.98-9.78) (3.29-4.09) (4.71-6.13) (0.00-2.58) 
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Figure 3 Median-joining network of mitochondrial cytochrome b haplotypes of Gekko gecko. The size of the circles is proportional to 


haplotype frequency. Each color represents locality. 


Table 5 Pairwise Fs; (lower diagonal) and P (above diagonal) 
values of Gekko gecko mitochondrial cytochrome b dataset. Bold 
values represent significant Fşr. 


Subclade Al A2 Bl B2 B3 
Al - 0.025 0.355 0.001 < 0.001 
A2 0.899 - 0.013 0.001 < 0.001 
Bl 0.989 0.971 - < 0.001 < 0.001 
B2 0.925 0.932 0.882 — < 0.001 
B3 0.897 0.906 0.755 0.849 — 


3.3 Population genetic structure and demography 
Pairwise genetic distances (% p-distance) between clades 
are presented in Table 4. The results showed high mean 
genetic divergences between G. reevesii and G. gecko 
(8.56—-10.59%). The mean p-distances between the 
major clades (Clades A and B) of G. gecko populations 
in Thailand were high (8.18—9.31%). The mean p- 
distances between populations within clade A were low 
(0.18-2.62%), whereas the mean p-distances between 
populations within clade B were higher than those within 
clade A (0.00-5.45%) and geographic structuring was 


Table 6 Results of analysis of molecular variance (AMOVA) across and within major geographic regions of G. gecko based on mitochondrial 


cytochrome b. 


Sources of Variation af. Sum of square Variance component Percentage of variation Fixation index P-value 
Among groups 1 658.71 25.28 46.24 0.46 0.099 
Among Pepe ys 3 1451.23 24.99 45.72 0.85 <0.001 
groups 

Within populations 109 478.87 4.39 8.04 0.92 < 0.001 
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detected, with the presence of three subclades (B1—B3). 

Overall mean pairwise Fr varied from 0.755 to 0.925 
and analyses revealed significant genetic differentiation 
between most subclades of G. gecko (Table 5). The high 
Fs, values indicated long-term restriction of gene flow 
between populations. AMOVA that revealed genetic 
differentiation among populations within group explained 
45.72% of the total genetic variation, and 8.04% of the 
genetic differentiation occurred within populations (Table 
6). 

The mismatch distribution was relatively ragged and 
multimodal in each subclade (Figure 4; Table 7). These 
tendencies were also found in Fu’s test of neutrality, 
showing non-significant values, whereas Tajima’s D 
was significantly negative (P = 0.011) only in subclade 
A2, rejecting the hypothesis of constant population size. 
However, this group contained individuals from different 
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regions in Thailand (CT, PL, NR, and SK). 
4. Discussion 


Whereas previous phylogenetic and phylogeographic 
studies of G. gecko have been conducted, they have 
included limited sampling from Thailand. This study 
carefully examined cyt b gene sequence diversity within 
and among 17 populations of G. gecko from residential 
areas in Thailand. Results of our analyses of partial 
cyt b gene sequences (1 125 bp) revealed two major 
mitochondrial clades (A and B) with high mean genetic 
divergence (8.18—9.31%). The mean genetic divergence 
within clade A ranged from 0.18 to 2.62% whereas that 
within clade B was higher, ranging from 0.00 to 5.45%, 
demonstrating that the genetic divergence among Thai 
populations in this study was much higher than those 


Table 7 Results of statistical test of neutrality and mismatch distribution of sampled populations of G. gecko based on mitochondrial 


cytochrome b. Bold values represent significant value. 


Clade Tajima’s D Fu’s Fs Mismatch Distribution 
D P Fs P SSD P r P 
Al’ - - -— -— - - -— -— 
A2 -1.6236 0.011 0.7516 0.644 0.6259 0.000 0.0907 1.000 
BI’ - - ~ — -— - -— -— 
B2 —0.08146 0.545 0.76685 0.691 0.01916 0.204 0.01921 0.544 
0.261 —0.50022 0.484 0.01611 0.164 0.02172 0.113 


B3 —0.75115 


Note: * means the subclade less than three individuals were excluded from the population demographic analyses. 


0.08 7 A Freq. Obs. 0.3 B —— Freq. Obs. 
T -—— Freq. Exp. =| -—— Freq. Exp. 
0.06 + J 
Jq 0.27 
0.04 : 1 
J 0.14 
0.02 AW a af 
ol er) ena A o 7 “~S a 
0 25 Wo o T7 100 125 0 10 15 20 
Pairwise differences Pairwise differences 
0.20 
E -€ Freq. Obs. Freq. Obs. 
a = Brog, Exp. -—— Freq. Exp. 
0.15 
0.10 | 
0.05 4 
>) ~ S FCO — a 
o 20 25 o 5 10 15 20 25 30 35 40 


Pairwise differences 


Pairwise differences 


Figure 4 Observed frequencies of pairwise nucleotide differences among mitochondrial sequences (dashed lines) and expected frequencies 
under a model of sudden population expansion (solid lines). Mismatch distributions depict frequencies of pairwise differences for: (A) All 
subclades of Gekko gecko, (B) subclade A2, (C) subclade B2 and (D) subclade B3. 
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previously reported between populations of G. gecko 
from China (Guangxi), Laos, and Vietnam (0.5—2.2% 
(424 bp) in Qin et al., 2007; 0.12-1.66% (1 141 bp) in 
Qin et al., 2012). High levels of mitochondrial genetic 
divergence were reported in other vertebrates in Thailand 
such as Hoplobatrachus rugulosus (cyt b; Pansook et 
al.2012), Chiromantis hasenae (16S rRNA; Yodthong et 
al., 2014), and Calotes mystaceus (COT, Saijuntha et al., 
2017), suggesting that strong genetic structuring of these 
species could resulted from geographic topography in 
Thailand and/or historical biogeographic process. 

The results of the median joining network analysis 
(Figure 4) were concordant with those of the phylogenetic 
analysis of G. gecko populations in Thailand, showing 
two major clades (A and B) and genetic sub-structuring 
within each clade. Clade B was divided into three 
subclades, B1 (ST population), B2 (Northeast, Central, 
and two individuals from RB2), and B3 (West, South, and 
one individual from CB). Strong genetic differentiation 
(Fsr) indicated limited gene flow among the subclades 
of geckos (Table 5) that could have resulted from the 
presence of geographic barriers, such as mountain ranges, 
preventing dispersal and gene flow between subclades B2 
and B3. The geographic topography of western Thailand 
contains steep hills and mountain ranges, i.e., Tenasserim 
and Thanon Thong Chai that run north-south orientation 
whereas the topography of northeastern Thailand 
contains plateau with lower elevations and mountain 
ranges, i.e., Phetchabun and Dong Paya Yen mountain 
ranges (Inger, 1999 ; Gupta, 2005). Therefore, these 
mountain ranges could play an important role as barriers 
of gene flow among G. gecko populations. A similar 
phylogeographic pattern was found in the rhacophorid, 
Chiromantis hansenae, with populations in northeastern 
Thailand separated from populations in western Thailand 
(Yodthong et al., 2014). These limited levels of gene flow 
could resulted from geographic barriers in western and 
northeastern Thailand and/or the limited dispersal ability 
of C. hansenae. In addition, the genetic structure in gecko 
populations could be affected by natural dispersal which 
is relatively limited in other geckos such as Hesperoedura 
reticulata and Hemidactylus mabouia (Hoehn et al., 2007; 
Short and Petren, 2011). Hoehn et al. (2007) suggested 
that a distance as small as 500 m is a barrier to the natural 
dispersal of H. reticulata. Short and Petren (2011) showed 
that genetic structure was found among populations of 
H. mabouia at a small regional scale and suggested that 
gene flow might be limited by both dispersal ability and 
geographic distance. Long-distance dispersal of some 
geckos e.g., H. mabouia, and Hemidactylus brookii has 


been reported in previous studies and could be explained 
by human-mediated transport (Short and Petren, 2011; 
Weterings and Vetter, 2018). Similarly, evidence of long- 
distance dispersal by human-mediated transport was 
revealed among populations of G. gecko in this study. For 
instance, H10 and H11 were shared by NR populations 
(northeastern Thailand), and RB2 population (western 
Thailand) whereas H5 was shared by CT (eastern 
Thailand), NR (northeastern Thailand), and PL (central 
Thailand) populations. Although the dispersal ability of 
G. gecko might affect the genetic structure of populations 
in Thailand, the dispersal pattern of this widely 
distributed species has not been studied. There is only a 
report on the foraging distance of G. gecko which ranges 
from 0.1—38.5 m from its retreat (Aowphol et al., 2006). 
Therefore, the dispersal patterns of G. gecko should be 
further investigated at the regional and landscape scales. 

Kongbuntad et al. (2016) analyzed the population 
genetic structure of G. gecko from northern and 
northeastern Thailand, Laos, and Cambodia based on 
multilocus enzyme electrophoresis. They found that G. 
gecko populations in northeastern Thailand were more 
closely related to populations from Laos and Cambodia 
and suggest that the populations from the Khorat Plateau 
of northeastern Thailand might be separated from 
northern populations by the Dong Paya Yen mountain 
range. In other lizard such as Calotes mystaceus, the 
phylogeographic patterns from most regions in Thailand, 
excluding the southern region, exhibited two lineages: (I) 
northern, western and central regions, (II) northeastern 
and eastern regions; and this divergence may be due to the 
Phetchabun, Dong Paya Yen, and Sankambeng mountain 
ranges acting as barriers to gene flow (Saijuntha et al., 
2017). These phylogeographic patterns are congruent 
with results from studies of other terrestrial organisms 
such as the centipede Scolopendra dehaani in mainland 
Southeast Asia (Siriwut et al., 2015). With these results 
and the findings of this study, sufficient evidence has 
been collected to suggest that intraspecific structure of 
terrestrial organisms in Thailand are primarily the results 
of historical biogeographic events in the region, namely 
the uplift of the Khorat Plateau (Racey, 2009). 

Although the sampled populations of G. gecko 
in Thailand have a shared evolutionary history, a 
significant amount of time would be required for 
different populations to establish unique haplotypes, 
suggesting independent evolutionary histories. The 
unique haplotypes indicate limited gene flow among 
populations within clades A and B, despite this species 
being well-adapted to the human-mediated environment 
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and the lack of obvious geographic barriers between 
some populations. Interestingly, the genetic divergence 
within clade A was low (0.18—2.62%). The subclade Al 
contained individuals from PL and RB2 whereas subclade 
A2 comprised populations from Laos, Vietnam, China, 
and eight individuals from NR, CT, PL, and SK. This 
finding indicates that G. gecko may have been introduced 
to different areas by anthropogenic activities such as 
commercial trade and human-mediated transport. In 
conservation assessments completed by the International 
Union for the Conservation of Nature (IUCN), G. gecko is 
listed as “least concern” because of its wide distribution, 
presumed large population and occurrence in natural 
and artificial environments; however, its population size 
and trends have not been evaluated (Lwin et al., 2019). 
International trafficking in geckos, especially G. gecko 
for both the pet trade and medicinal uses, occurs on a 
grand scale (Bauer, 2009). Lwin et al. (2019) suggested 
that international trade monitoring is necessary, including 
CITES monitoring to collect data on trade volumes. 
In conclusion, the results obtained through this study 
demonstrated substantial genetic differentiation among the 
G. gecko populations in Thailand, indicating limited gene 
flow that could be caused by the complexity of Thailand’s 
topography and the dispersal ability of G. gecko. We 
recommend the extension of fine-scale sampling across 
the species’ distribution ranges, combining morphological, 
and molecular investigations in future studies to facilitate 
understanding of the evolutionary history, and the 
demographic and taxonomic status of G. gecko (i.e., G. g. 
gecko and G. g. azhari). Additionally, this study provided 
evidence that human-mediated transport has affected the 
population genetic structure of G. gecko populations in 
mainland Southeast Asia (Thailand, Laos and Vietnam) 
and southern China. Therefore, conservation strategies 
should involve controlling the introduction pathways 
of non-native G. gecko populations such as commercial 
and pet trade routes to prevent further impacts on native 
populations and other species. 
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